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ABSTRACT 

The surface brightness diameter relationship for supernovae remnants (SNRs) 
is explained by adopting a model of direct conversion of the flux of kinetic en- 
ergy into synchrotron luminosity. Two laws of motion are adopted, a power law 
model for the radius-time relationship, and a model which uses the thin layer 
approximation. The fluctuations on the log-log surface diameter relationship are 
modeled by a Monte Carlo simulation. In this model a new probability density 
function for the density as function of the galactic height is introduced. 

Subject headings: ISM: supernova remnants, radio continuum: galaxies 



Introduction 



The correlation between radio surface brightness, £iae>> an d diameter D in pc in su- 
pernova remnants is parametrized as 
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wher e /3 E and Cs are are found through the radio observations (I S hklovskiil 1 1960 : lUrosevic et al 
20051 ) . This is an observa t ional relationship which has been explained by early theoreti 
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1979 ). A theoretica l interpretation of the £ — D relationship has been introduced 



by iDuric &: Seaquistl ( 119861 ) by incorporating the Sedov blast wave solution for supernovae 
remnant (S NR) expansion, the generation and evolution of the magnetic field as outline d by 
and the acceleration of relativistic electrons by shocks, as formulated by iBelll 
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(1978a 
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bj). The time-dependent, nonlinear kinetic theory for cosmic ray (CR) production was 
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adopted by iBerezhko fc York! (120041 ) in order to derive a theoreti cally predicted brightness- 
diameter relation in the radio range for the Sedov phase. Recently iBandiera fc Petrukl (120101 ) 
found that SNRs cease to emit effectively in radio at a stage near the end of their Sedov 
evolution, and that models of synchrotron emission with constant efficiencies in particle ac- 
celeration and magnetic field amplification do not provide a close match to the data. These 
previous models leave a series of questions unanswered or merely partially answered: 



• Is it possible to deduce a theoretical £ — D relationship on the basis of power law 
expansion for SNR? 

• Is it possible to deduce a classical equation of expansion for the SNR with two ad- 
justable parameters that can be found from a numerical analysis of the radius-time 
relationship? 

• Is the theoretical H — D relationship sensitive to the density of the interstellar medium, 
which decreases as a function of height above the galactic plane? 



In order to answer these questions, we develop a theoretical derivation of the S— D relationship 
by modeling the expansion of a synchrotron emitting shell. In Section (2J we show that it is 
reasonable to fit a power-law expansion relationship to a known young supernova remnant, 
SN 1987J in M81, for which the growth of the shell like structure is available over a period 
of 10 years. In Section [3j we discuss a modification of the power-law expansion model for 
the effects of the standard gradient in density in the interstellar medium. In Section HI we 
complete our derivation of the E — D relationship by adopting a direct conversion of the 
flux of kinetic energy into radiation calibrated on the observational data (Section ED- In 
Section (6J we address the observed fluctuations in the £ — D relationship by utilizing a new 
probability density function (PDF) for the contrast of in situ densities of SNRs. This new 
PDF is derived from the standard number distribution of SNRs with the galactic height. 



2. The Power Law Model of Expansion 

A standard point of view assumes that the SNR expands at a constant velocity until 
the surrounding mass is of the order of the solar mass. This time scale can be modeled by 

t M = 186.45 - yr , (2) 

\Pk>v wooo 

where M is the number of solar masses in the volume occupied by the SNR, n , the number 
density expressed in particles cm -3 , and fioooo the initial velocity expressed in units of 



-3 - 



10000 km/s , see iMcCrayl ( 119871 ). To model the time- dependent expansion of SNRs, we use 

R(t) = Ct a , (3) 

where the two parameter s C and a can be f ound from the observations. Ten years of 
observations of SN 1993 J (IMarcaide et al.ll2009l ) allow us to fit these two parameters which 
are reported in Table [T] and in Figure [TJ The quality of the fit is measured by the merit 
function 
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where Rth, Robs and a \, s are the theoretical radius, the observed radius and the observed 
uncertainty respectively. This observed relationship allows us to express the radius, as a 
function of time, by 

R(t) = 0.0155 1 - 828 pc , (5) 
where the time t is expressed in yr. The velocity in this model is 

km 



V(t) = 12587.67r am 



(6) 



Figure [2] displays the observed instantaneous velocity as deduced from the finite difference 
method and the best fit according to Eq. ([6]). Another interesting SNR 
started to be visible in 1006 AD and has an actual radius of 15' according to 



The r adius in pc is a function of the the adopted distance; a s an example, iKatsuda et al. 
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Katsuda et al. 



(120101 ) quotes a distance of 2.2 kpc and IWinkler et al.l (120031 ) 2.18 kpc. Here we adopt a 
distance of 2.2 kpc and consequently the observed radius, R b s , is 



R = 9.59D 22 pc 



(7) 



where the D22 is the distance expressed in units of 2.2 kpc. The proper motion of the 
remnant's edge is a weak function of the azimuth angle and the average value is 



Fobs = 0.5 arcsec/yr 
The average velocity of the expansion is 

Vobs 
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Eq. ([3]) and the connected velocity can be solved for the two unknown variables a and C 

V bs t 
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where V oos is expressed in — , R obs in pc and t in yr. Table [2] reports the numerical results. 
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3. Momentum conservation and galactic height 

The presence of a slow wind from the SN progenitor makes it reasonable for us to assume 
that the SNR evolves in a previously ejected medium phase in which density is considerably 
higher than the interstellar medium (ISM). The assumption used here is that the density of 
the ISM around the SNR has the following two piecewise dependencies 

where / is a parameter that models the jump in density , < / < 1, and is connected with 
the vertical profile in density as a function of the galactic height. In this framework, the 
density decreases as an inverse power law with an exponent d that can be fixed from the 
observed temporal evolution of the radius, with d = meaning constant density; further on 
the parameter / regulates the density at R — Ro- The swept mass in the interval < r < i? 

iS 4 

M Q = -p 7iR 3 . (12) 

The swept mass in the interval < r < R with r > R is 

M 
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PottRq 3 4 3 
+4/ d _ 3 +-p«kR g . (13) 

Momentum conservation requires that 

Mv = M v , (14) 

where v is the velocity at t and v is the velocity at t — t . The velocity as a function of the 
radius is 

3r d R3~ d f + r 3 (-d + 3-3f) ' 1 ' 

In this first order differential equation in R, the variables can be separated and an integration 
term by term yields the following nonlinear equation 

J~nl — 

(12 R 3 f - 3 R 3 fd + 7 R 3 d - R 3 d 2 - 12 R 3 ) R 
-3 R d fR^ d + 7R 3 v dt + R 3 v dH 
-R 3 v d% + R ^d 2 + 12 - 9 R 4 f - 7 R 4 d - 7 R 3 v dt 

-12 R 3 v t + 3 R 4 fd+ 12 R 3 v t . (16) 
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An approximate solution of Fnl(t) can be obtained assuming that 
3R$R 4 - d f > R 3 (4 - d) (d - 3 + 3 /) R 

R(t) = 

(Ro 4 ~ d -jj(d-3 + 3f)R 4 ~ d (A - d) 
+jj(A-d)v R 3 - d (3-d)(t-t ))^ . 



(17) 



The physical units have not been specified, pc for length and yr for time are perhaps 
an acceptable astrophysical choice. With these units, the initial velocity vq is expressed in 
p and should be converted into ^p; this means that v = 1.02 10~ 6 i;i where V\ is the initial 
velocity expressed in — . 

A one - dimensional solu tion of Eq. ( f!6l) can be found with the FORTRAN subroutine 
ZRIDDR (jPress et al.lll992l ). A plot of this radial solution, using a parameter value of f = 
1, is shown in Figure ([3]). Plots of other solutions with different value of f can be found in 
Figure (Hj). 



4. E — D in the power law model 

In this section, we review the time scale of synchrotron losses, the time scale of accel- 
eration and the inequality which allows us to use the in situ acceleration of electrons. The 
relation between the radio surface brightness and diameter of SNRs is deduced assuming a 
direct proportionality between the radio-luminosity and the flux of kinetic energy. 



4.1. The synchrotron emission 



An electron which loses its energy due to synchrotron radiation has a lifetime of 

E 

Tr ~ ~ bOOE^H^sec , (18) 
where E is the en ergy in ergs, H the magnetic field in Gauss, and P r is the to tal radiated 



power (]Langill999l . Eq. 1.157). The energy is connected to the critical frequency (]Langill999 
Eq. 1.154), as 

v c = 6.266 x 10 18 HE 2 Hz . (19) 
The lifetime for synchrotron losses is 

1 



r syn = 39660 
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yr 



(20) 
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Following ( jFermil Il949l . 119541 ). the gain in energy in a continuous form for a particle 
which spirals around a line of force is proportional to its energy, E, 



dE 
~dt 



E 

m 



where tjj is the typical time-scale, 
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Tjj 3 (T L 



(21) 



(22) 



where u is the velocity of the ac celerating c loud, c is the velocity of the light and L is the 
mean free path between clouds (jLangill999L Eq. 4.439). The mean free path between the 
accelerating clouds in the Fermi II mechanism can be found from the following inequality in 
time : 

tjj < T sync , (23) 
which corresponds to the following inequality for the mean free path between scatterers 



T 1.72 10 5 w 2 
L < pc 



HVHvc 2 

The mean free path length for SN 1993J at v — 1GHz gives 

L < 1.831(T 5 r a3 V 



(24) 



(25) 



where the velocity is given by Eq.()6]) and H = 65.1 G (IMartf-Vidal et al.ll201ll ). When this 
inequality is verified the direct conversion of the flux of kinetic energy into radiation can 
be adopted. Recall that the Fermi II mechanism produces an inverse power law spectrum 
in the energy of the type N(E) oc E~J or an inverse p ower law in the observed frequencies 
N(i>) oc v~ a with a 
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4.2. The dimensional approach 



The source of synchrotron luminosity is here assumed to be the flux of kinetic energy 

1 



L m = -pAV 3 



(26) 



where A is the emitting surface area (jde Young! 120021 . Eq. A28). Assuming the surface area 
is like that of a sphere, the luminosity is 



^ P AttR 2 V 3 



(27) 
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where R is the instantaneous radius of the SNR and p is the density in the advancing layer 
in which the synchrotron emission takes place. Assuming the density of the advancing layer 
scales as R~ d , which means that 

L m oc R 2 ~ d V 3 . (28) 
The time dependence is eliminated by utilizing Eq. fl3]) 

L m = L (—) — , (29) 

n 

where L is the luminosity at R = R® or 

_ T ., ID . da — 5 a + 3 , . 

L m = L (— ) — , (30) 

where D is the actual diameter and Dq = 2Rq. On assuming that the observed luminosity 
,L U , in a given band denoted by the frequency v is proportional to the mechanical luminosity 
we obtain 

. D . da — 5 a + 3 t 9 / \ 

Lu = L ^{—) - Jykpc , (31) 

where Lq iU is the observed radio luminosity in a given band. The observations are generally 
represented in the form 

L = C L D-to , (32) 



where L and Cj, are parameters deduced from the radio observations (IGuseinov et al.l 12003 



Urosevic et al.l 120051 1 . Given the two observational parameters a and /3l we can derive 



d= ^ + 5 a -3 . (33) 
a 



The radio surface brightness of a remnant is defined as 

v _ SiGHz 

^-—ffT- > 



(34) 



where S\guz is the detected flux of a remnant at 1GHz and 9 the observed angle. Due to 
the fact that 9 oc D 2 we have 

^G Hz = ^, lGHz {—) - , (35) 

where £ ,igh2 * s ^ ne sur f ace brightness at -D = D . According to the basic observational 
relationship (JTJ 

^= fea + 3Q - 3 . (36) 
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5. Observations in the power law model 



Supernova remnants in our galaxy have two luminosity-diameter relationships (jGuseinov et al 

20031 ): the first one for SNRs which have L > 5300 Jy kpc 2 , D < 36.5 pc and the second one 



for SNRs having L < 5300 Jy kpc 2 , D > 36.5 pc, 

L = 2.45 10 4 L>-°' 43 Jy kpc 2 D < 36.5pc , 

and 

L = 5.38 10 9 /r 3 - 84 Jy kpc 2 D > 36.5pc . 
The E — D relationship in a similar way is 

E = 2.7+ 2 ;* 10- 17 D~ 2A7 -o ™ W mT 2 R z~ x stcr' 1 D < 36.5pc 

and 



8.4^ 9 3 5 10- 12 D- 5 - 99 -°£ Wm- 2 Hz- l stcr- 1 



D > 36.5pc 



(37) 
(38) 

(39) 
(40) 



The observations of t he surface brightness of SNRs in other galaxies has been analyzed by 
Urosevic et al.l (120051 ) with data available at Centre de Donnes astronomiques de Strasbourg 
(CDS). Figure |5] reports the observational data as well the fitting curve of all the radio SNRs 
in external galaxies. The results of our numerical analysis for extragalactic SNRs gives 



E = (8.8 ± 3.27 Kr 16 zr 3 - 1±0 - n Wm~ 2 H zoster' 1 D < A50pc 



(41) 



which agrees with the results of lUrosevic et al.l (120051 ) . The values of d from the 
E — D relationship are reported in Table |3j 



6. A Monte Carlo model for E — D 

The significant fluctuations observed in the E — D relationship can be explained by 
correlating the well known probability to have a SNR at the galactic height z with the 
density. This conversion can be done using the nonlinear relationship between z and density 
as given by the self-gravitating disk. 

The radius function we derived in Section |3] (Eq. H71) from momentum conservation 
allows us to deduce an analytical expression for the E — D relationship: 

Eicii* = ^o,ighz(-^)(jz-) 5 • ( 42 ) 

j 6 Dq m z nzsr 



The E — D relationship now has an inverse cubic dependence for the contrast parameter /. 
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6.1. The profile of the ISM 



The vertical number density distribution of galactic H I has the following three compo- 
nent behavior as a function of the galactic height z , the distance from the galactic plane in 
pc : 

n(z) = n ie - z2/Hl2 + n 2 e~ z2/H22 + n 3 e-^ /H3 . (43) 

We set the densities in Eq. (1431) to n\ = 0.395 particles cm -3 , n 2 = 0.107 particles cm -3 , 
n 3 = 0.06 4 particles cm" 3 , and the scale heights to = 127 pc, H 2 = 318 p c, and 
# 3 = 403 jLockmanill98i ; bickev fc Lockman![l9~9ol : iBisnovatvi-Kogan fc Silich] Il995h . This 



distribution of galactic HI is valid in the range 0.4 < R < Rq, where Rq = 8.5 kpc and R 
is the distance from the galaxy center. The previous empirical relationship can be modelled 
by a theoretical one. The density profile of a thin self-gravitating disk of gas which is 
characterized by a Maxwellian distribution in velocity and distribution which varies only in 
the z-direction has the following number density distribution 



n(z) = nosech 2 



' 2 * z ' 



(44) 



wher e no is the density at z = and z Q is a scaling parameter (IBertinl 120001 ; iPadmanabhan 
20021 ). Figure ([6]) displays a comparison between the empirical function sum of three expo- 
nential disks and the theoretical function as given by the Eq. ( 1441) . The previous equation 
can be expressed with our contrast parameter / 



/ = sech 2 



and the inversion gives z as function of / 



1 2 * z ' 



2 arcsech y\fj 



(45) 



(46) 



6.2. The distribution of SNRs 



The probability density function (PDF) ,p(z), to have a SNR as function of the galactic 
height z is characterized by an exponential PDF 



p(z) 



with b = 83 pc ( jXu et al.l 120051 ). The average value z for SNRs is 

z = b = 83pc . 



(47) 



(4* 



- 10 - 



We briefly review how a PDF p(z) changes to g(f) when a new variable f(z) is introduced. 
The rule for transforming a PDF is 

^)- p m • (49) 

dz 

Once the previous rule is implemented we obtain 

_2 arcsech ( -/J) z 

9(f) = ~^S= . (50) 



V 3/ V^ - V77 + 1 
The average value of / when b = 83 pc and Zq=90 is / = 0.7927. 



6.3. Monte Carlo Simulation 

We are ready to build a a Monte Carlo model for the S — D distribution of SNRs for 
SNRs characterized by the following constraints 

• The lifetime of a SNR is generated between and t max . 

• The time is converted in diameter , D = R * 2 , adopting Eq. ([5]). 

• A value of / is randomly generated according to the PDF, Eq. (l50]) . 

• In the previous points we have generated one D and one /, as a consequence a value 
of T,i GHz as given by Eq. (H2|) can be generated. 

The results of this simulation are displayed in Figure [7] for the Galactic distribution, and in 
Figure [S]for the extragalactic distribution. In this Section we derived the profile of the density 
in our galaxy as function of the galactic height z in the framework of the self-gravitating 
disk. This new relationship allows us to find a new PDF that characterizes the effect of 
SNRs expanding into a medium of varying density. The great fluctuations present in the 
galactic S — D relationship are therefore explained by the SNR evolution in a medium with 
a lower density. 



7. Conclusions 



The relation between radio surface brightness and diameter of galactic and extragalactic 
SNRs presents a linear relationship in the Log( S)-D plane. Superposed on this observed 
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linear behavior there is a fluctuation which is due to the variable density with the galactic 
height of the ISM. The major results of this paper are summarized as follows. 

• The power law model for the expansion of a SNR coupled with a decreasing density 
medium allows to build a simple expression for the luminosity as well the radio surface 
brightness, see Eqs. (135]) and ( 129]) . 

• A more sophisticated approach to the law of motion as given by the momentum con- 
servation in a medium with variable density allows to obtain the same result , see Eq. 
( )42l) . once an asymptotic behavior for the law of motion is given, see Eq. ( fl71) . In this 
model is possible to introduce a contrast parameter / which allows to state that both 
the radius and the velocity increase as 4 and the surface brightness as -p. 

• The value of / for each SNR is given by a new PDF (Eq. [501) which is derived from 
the observed PDF for SNR as function of the galactic height , see Eq. (j4"T|) and from 
a physical profile of density of the ISM as given by Eq. (1441) . 

• A Monte Carlo simulation of the radio surface brightness explains the fluctuations 
visible in the X — D relationship as a probabilistic effect to have a SNR with a given 
galactic height or /. 
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Fig. 1. — Theoretical radius as given by the power law model with data as in Table [T] (full 
line), and astronomical data of SN1993J with vertical error bars. 




Fig. 2. — Theoretical velocity as given by Eq. ([6]) (full line) and instantaneous velocity of 
SN 1993J with uncertainty. 
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Fig. 3. — Theoretical radius as obtained by the solution of the nonlinear Eq. (TlBT) (full 
line), data as in Tabled] and f=l. The astronomical data of SN 1993J are represented with 
vertical error bars. 
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2 4 6 8 10 

Time (yr) 

Fig. 4. — The same as Figure [3] but / is variable; / = 1 (full line), / = 0.9 (dashed) , 
/ = 0.8 (dot-dash-dot-dash) / = 0.7 (dotted) and experimental point represented by empty 
stars. 



-18 - 




Fig. 5. — The L OG-LOG S — D diagr am at a frequency of 1 GHz for all the extragalactic 
SNRs as given by lUrosevic et all (120051 ). 
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Fig. 6. — Profiles of density versus scale height z: the self-gravitating disk as given by Eq. 
(j44p when zq = 90 pc (dashed) and the three-component exponential distribution as given 
by Eq. (g3D (full line). 
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log10(D(pc)) 



Fig. 7. — A Monte Carlo LOG-LOG Y, — D diagram at a frequency of 1 GHz which simulates 
the Galactic SNRs for D < 36.5pc. The parameters are t max — 5055 yr , a = 0.82 , b = 83 
pc and z = 90 pc. The Monte Carlo simulation gives C s = 2.7110 -17 and /3 S =2.45. 
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Fig. 8. — A Monte Carlo LOG-LO G E — D diagram at a frequency of 1 GHz (green circles) 
and extragalactic SNRs as given by lUrosevic et al.1 ( 120051 ) for D < 36.5pc (empty red stars). 
The parameters are t max = 2960yr , C = 0.17 , a = 0.57 , b = 83 pc and zo = 60 pc. The 
Monte Carlo simulation gives = 8.8410 -16 and /?s=3.65. 
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Table 1: Numerical values of the parameters of the fits for SN 1993J and x 2 - N represents 
the number of free parameters. 



N values \ 2 

power law 
2 a = 0.82 ± 0.0048 6364 

C = (0.015 ± 0.00011) 
nonlinear radius 
4 d=2.93; r = 0.019 pc; 276 
t =0.249 yr; v =100 000 — 
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Table 2: Numerical values of the parameters of the fits for SN1006 as given by Eq. (fTUj) 
for R obs = 9.59 pc , V obs =5441 *p and t= 998 yr. 



unknown 


value 


a 


0.579 


C 


0.1758 
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Table 3: Numerical values of d 



environment 


Pobs 


a=0.82 


a=0.57 


extragalactic, D < 450 pc 


3.1 


2.44 


0.83 


galactic , D < 36.5 pc 


2.47 


1.81 


0.2 



